Elastic Theory Has Zero Radius of Convergence 



Alex Buchel and James P. Sethna 
Laboratory of Atomic and Solid State Physics, 
Cornell University, Ithaca, NY, 14853-2501 



< 



o 
o 



Nonlinear elastic theory studies the elastic constants of a material (such as Young's modulus or 
bulk modulus) as a power series in the applied load. The inverse bulk modulus K, for example 
depends on the compression P: 1/K(P) = Co + ciP + C2P 2 ■ ■ ■ + c n P n + ■ ■ ■. Elastic materials 
that allow cracks are unstable at finite temperature with respect to fracture under a stretching 
load; as a result, the above power series has zero radius of convergence and thus can at best be an 
asymptotic series. Considering thermal nucleation of cracks in a two-dimensional isotropic, linear- 
elastic material at finite temperature we compute the asymptotic form c n +i/c n — > Cn ' as n — > 00. 
We present an explicit formula for C as a function of temperature and material properties. 
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Hooke's law, F = kx, representing the elastic response of a body to an external stress, is only the first term in a 
Taylor's series. This paper addresses the convergence properties of this expansion. We argue generally that the radius 
OO , of convergence of elastic theory is zero, and for two-dimensional isotropic linear-elastic theory allowing for brittle 
t-H ' fracture, we can calculate the asymptotic behavior of the coefficients. 
For simplicity, let's consider the bulk modulus K(P), 

Under stretching (P < 0), the true ground state is fractured into pieces (relieving the strain energy). As a result, 
f^*) ■ P = can not be a point of analyticity for K(P), and thus (1) has zero radius of convergence. 

Similar arguments were used by Dyson in 1952, where he argued that calculations in quantum electrody- 
namics, expressed as a power series in the fine structure constant e 2 /frc ~ 1/137, had zero radius of conver- 
gence (because negative values of e 2 lead to unstable theories). This did not prevent these calculations from be- 
ing useful (indeed, they represent the best quantitative agreement between theory and experiment known to sci- 
ence). The community believes these expansions are asymptotic in the same sense as Stirling's approximation 

T(n) = (n - 1)! ~ e~ n n n ^27r/nj ^1 + l/(12n) + l/(288n 2 ) + • • : at any fixed n no matter how large, Stir- 
ling's series in 1/n eventually diverges, but the difference between the function and the M th approximation goes to 
zero faster than l/n M as n — > 00. 
, I Since Dyson's work, field theoretic methods have been developed to relate the instabilities in the theories at 
small negative couplings to the high-order terms in perturbation theory. Here we apply these methods to a particular 
case, using the statistical mechanics of thermally nucleated cracks to calculate the high-order terms c n in the inverse 
bulk modulus (1). Statistical and thermodynamic approaches to crack nucleation and fracture have an established 
history However, most work in this area is concentrated on failure at rather high stresses, near the threshold for 
instability (the spinodal point). The high-order terms in the perturbative expansion of the inverse bulk modulus are 
governed by the elastic response of the material to infinitesimally small tension (see below) , so we are far away from 
the spinodal point and linear elastic theory is an adequate description. 

Consider an infinite two-dimensional isotropic linear-elastic material subject to uniform compression P at infinity. 
Creation of a cut of length i will increase the energy by 2at, where a is the surface tension (the energy per unit length 
of edge), with a factor of 2 because of the two free surfaces. On the other hand, for negative P (uniform tension) the 
cut will open up because of elastic relaxation. Calculating this relaxation energy we find the total energy E of a crack 
of length I: 

■kP 2 (1 - <7 2 )£ 2 

E(£)=2ae- ^- K \ Y a)i . (2) 

Introducing 
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4 = nP 2 {l~a 2 ) (3) 
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we can rewrite the energy of the crack as 



£ 2 

E = 2a£-a—. (4) 

It follows that cracks with £ > i c will grow, giving rise to the fracture of the material, while those with £ < £ c 
will heal — a result first obtained by Griffith |Q. This is the instability that is responsible for the breakdown of 
elastic perturbation theory. Because the energy E{£ c ) = a£ c grows as 1/P 2 as P — > 0, interactions between thermally 
nucleated cracks are unimportant at small P and low temperatures (allowing us to use the "dilute gas approximation"). 

The thermodynamic properties of a macroscopic system can be obtained from its partition function Z: 

Z = J2^M-/3E n ) (5) 

n 

where the sum is over all states of the system. Once a perturbative expansion for the free energy F = hiZ is 

known, one can calculate the power series expansion for the inverse bulk modulus using 



1 1 /OF 

K(P) ~ ~PA\dP 



(6) 



where the elastic material has area A. 

For P < 0, our model is in a metastable state, and direct computation of the partition function should yield 
a divergent result. A similar problem for the three-dimensional Ising model was solved by Langer ||: one has to 
compute the partition function in a stable state P > 0, and then do an analytical continuation in parameter space to 
the state of interest. The free energy develops an imaginary part in the unstable state, related to the decay rate for 
fracture |J : the situation is similar to that of barrier tunneling in quantum mechanics ]To[ | , where the imaginary part 
in the energy gives the decay rate of a resonance. 

The calculation of the imaginary part of the partition function is dominated by a saddle point, that in our case is 
a straight cut of length £ c . The straight cut is the saddle point because it gains the most elastic relaxation energy for 
a given number of broken bonds (as can be checked with a direct calculation [pT[). 

For simplicity, we start by considering the model without including the quadratic fluctuations around the saddle 
point. The partition function for a dilute gas of straight cuts of arbitrary length at their equilibrium shape for the 
tension P can be calculated directly: the imaginary part of the free energy for P < is 

--•<-)^(W^y) V2 K)-{^)} 

with A being the ultraviolet cutoff in the theory (roughly, the interatomic distance). (The factor (2ttA/X 2 ) comes 
from the zero modes for rotating and translating the cracks.) 

Assuming the free energy is analytic in the complex P plane except for a branch cut along P G (— oo, 0], we obtain 
a Cauchy representation for the free energy 



As was first established for similar problems in field theory |2|-||, this relation determines the high-order terms in 
perturbative expansion of the free energy F(P) — J2 n fn,P n 



7T J P' n+1 

— oo 

Because the saddle point calculation becomes more and more accurate as P — ► 0, and because the integrals in 
equation (9) are dominated by small P as n — ► oo, using the saddle-point form for the imaginary part of the free 
energy yields the correct n — > oo asymptotic behavior of the high-order coefficients /„ in the free energy. For our 
initial approximation (7) 
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/rw = ( _ 1) ,„ lr ^_^ i __ i ^ ^Jss^s do) 

We can then use the thermodynamic relation (6) to show c„ = — (n + 2) f n + 2 / 'A, and thus calculate the asymptotic 
behavior of the expansion of the bulk modulus: 



Cn+l _^ _„l/2( ^i 1 - (j2 ) 
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" 1/2 ( oa Z" 2 ) asn->oo (11) 



which indicates that the high-order terms c n in the perturbation expansion for the inverse bulk modulus roughly 
grow as (n/2)!. The formula (11) will still be correct once we add the quadratic fluctuations (although a will develop 
temperature-dependent corrections). For cracks in three dimensions, we expect the scaling c n +i/c n ~ n 1 / 4 as n — > oo, 
using similar arguments. 

The above calculation ignores the quadratic fluctuations about the saddle point, which will change the prefactor in 
the expression (7) for the imaginary part of the free energy and will renormalize the surface tension a. There are two 
kinds of quadratic fluctuations we have to deal with. (I) Curvy cuts — changes in the shape of the tear in the material: 
deviations of the broken bonds from a straight-line configuration. (II) Surface phonons — thermal fluctuations of the 
free surface of the crack about its equilibrium opening shape. Since the energy of the curvy cracks is calculated for 
their equilibrium configuration, the response of the surface phonons to curving the cut is already incorporated, so the 
quadratic fluctuations (II) can be calculated independently from (I). 

In both cases the theory needs regularization: we must decide exactly how to introduce the ultraviolet cut-off 
A. We describe the curvy crack as a bunch of line segments, parameterized by kink angles a, — angles between 
consecutive segments. We calculate the energy release in the material with a crack using a conformal mapping from 
the plane with a cut to the exterior of a unit circle [fijl, which for a piecewise linear cut is given by a modified 
Schwartz-Christoffel map. It is convenient for the regularization of the surface phonons to use equally spaced points 
on the unit circle, of angular spacing tt\/£; this corresponds to points which accumulate at the edges of the original 
crack cut. For consistency, we use the same regularization for the curvy crack segments. We have also used the (less 
convenient, more natural) regularization of equally spaced points along the original crack cut, and the results are 
similar (although the convergence is not as convincing): the form of the imaginary part of the free energy and the 
ratio of the high-order elastic coefficients (11) do not change, but (see below) the constants sq, si and s 2 and the 
temperature-dependence of a do change JO) . 

To calculate the quadratic fluctuations due to the surface phonons, we must calculate the determinant of the matrix 
M p , given by the energy spectrum of the normal modes about the crack's equlibrium (saddle point) opening: 

dctM > = 2 * (w 1 ) 3 - -» {" ( 2 k - 2 ) } < 12 > 

where n is number of kinks. Similarly, the quadratic fluctuations due to curving the crack is given by the determinant 
of the matrix My that specifies the change in energy release to quadratic order in the kink angles a.i (so the energy 
change due to the curving of a straight cut AE({cti}) = MijOtiCtj +0(a 3 )). We have obtained an analytical expression 
for the entries of the above matrix by calculating, to quadratic order, the energy of a crack with two kinks. This 
lengthly and tedious calculation |ll[] was checked both using a finite-element crack simulation |l3| and an exact 
analytical solution for a cut that is an arc of a circle. We assume that scaling of the curvy crack determinant is of the 
same form as for the surface phonons : 

detM ~ s n Sl exp {n\ln + s 2 ^j j (13) 

as the number of kinks n — > oo. We extract the constants so, Si, and s 2 numerically, using systems with up to n = 400 
kinks: we extrapolate visually to find lnso = 1000 ± 1000, s\ = —400 ± 100, and s 2 = 7.5 ± 0.5. Notice that sq is 
subdominant to si, which is subdominant to s 2l so our errors get rather large - even without including our uncertainty 
as to the scaling form (13). 

These quadratic fluctuations modify the imaginary part of the free energy (7) as follows: 



/8^1 - cr 2 V A 



-Kl( 1 -T + ta W + 3 h ^)} to ^ m 
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with l c given by (3) and lmF simple (P) given by (7). The first term in (14) modifies the prefactor, the other one 
effectively renormalizes the surface tension: 



1 f s 2 3YX 2 1 BaX\ 



Finally, we derive the correct asymptotic form for c 

t(1 - a 

A8Ya 2 r J \X 2 J\ a J 8\/27r/?a(l - a 2 ) 

This is the main result of our paper. 

We warn the reader to treat our results in the proper context. First, we do not expect this calculation to have 
experimental implications in the near future. Real fracture invariably occurs on inhomogeneities in the material: 
pre-existing surface or bulk microcracks, dislocation tangles, or grain boundaries. Even for a perfect dislocation-free 
crystal with stabilized surfaces, the effects we describe will remain immensely small: the reason that incredibly tiny 
compression make (relatively) large changes in volume (leading to large high-order terms in the power series) is because 
they suppress the (already incredibly rare) opening of large thermally nucleated cracks: measuring these effects would 
seem infeasible. 

Second, our results can be viewed as a straightforward extension to the solid-gas sublimation point of Langer [^]|| 
and Fisher's fli)] theory of the essential singularities at the liquid-gas transition. Indeed, if we allow for vapor pressure 
in our model, then our system will be in the gas phase at P = 0. The essential singularity we calculate shifts from 
P = to the vapor pressure. If we measure the nonlinear bulk modulus as an expansion about (say) atmospheric 
pressure, it could converge: but the radius of convergence would be bounded by the difference between the point of 
expansion and the vapor pressure. 

Third, we have forbidden dislocation nucleation and plastic flow in our model. Dislocation emission is crucial for 
ductile fracture, but by restricting ourselves to a brittle fracture of defect-free materials we have escaped many (lethal) 
complications. Dislocations are in principle important: the nucleation |l5f| barrier Edis for two edge dislocations in 
an isotropic linear-elastic material under uniform tension P with equal and opposite Burger's vectors b is 



Yb 2 , Y 
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4tt(1-o- 2 ) P 



E ^s = — 2T ln o + E o (17) 



where Eq is a P independent part that includes the dislocation core energy. The fact that Edis grows like 1/ ln P as 
P — > (much more slowly than the corresponding barrier for cracks) tells that in more realistic models dislocations 
and the resulting plastic flow (ITJ] cannot be ignored. While dislocations may not themselves lead to a catastrophic 
instability in the theory (and thus to an imaginary part in the free energy?), they will strongly affect the dynamics 
of crack nucleation as noted above. 

Fourth, we are being careful to distinguish tension from negative hydrostatic pressure. The response of a crystal to 
a force depends in a fundamental way on whether the force couples to the atoms or to the lattice (i.e., to the broken 
symmetry). Under gravity, all crystals will flow like anisotropic liquids, with a viscosity given simply in terms of the 
vacancy diffusion constant. In contrast, a crystal gripped at the sides and strained will not flow with a rate linear 
in the external force: under a strain produced by coupling to the surface layers, single vacancies cannot individually 
relieve the strain, and a dislocation pair (a loop in a three-dimensional case) must appear in the material for strain 
relief to occur: the plastic response of a two dimensional crystal to external strain is a temperature-dependent power- 
law in the external force |fl6|| . Although we call the compression (or tension) P, it is not hydrostatic pressure, but 
a coupling to the lattice. Under negative hydrostatic pressure, vacancies will have a negative chemical potential /i, 
and the dominant fracture mechanism becomes the nucleation of vacancy clusters or voids (rather than Griffith-type 
critical microcracks), as noted by Golubovic and collaborators [r^| . For negative chemical potential /i per unit area, 
their vacancy clusters (up to constants) have energy E vac (R) ~ aR — |^|i? 2 . If we identify \x with P, comparing with 
(2) we see that the vacancy cluster gains an energy linear in P, while the crack gain is quadratic. This leads to a 
problem that maps onto Langer's calculation, leading to the asymptotic relation c n+ i/c n ~ n, so the coefficients in 
this case would diverge more strongly: c n ~ nl. However, the identification of [i with P demands a mechanism for 
relieving elastic tension by the creation of vacancies. Ignoring possible effects of the boundaries (presumed infinitely 
far removed), as noted above vacancies will be created (and hence relieve tension) only through dislocation motion, 
which we have explicitly excluded from our model. 
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